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Abstract. We consider an equation for a thin-film of fluid on a rotating 
cylinder and present several new analytical and numerical results on steady 
state solutions. First, we provide an elementary proof that both weak and 
classical steady states must be strictly positive so long as the speed of rotation 
is nonzero. Next, we formulate an iterative spectral algorithm for computing 
these steady states. Finally, we explore a non-existence inequality for steady 
state solutions from the recent work of Chugunova, Pugh, &; Taranets. 



1. Introduction 

Thin liquid films appear in a wide range of natural and industrial applications, 
such as coating processes, and are generically characterized by a small ratio between 
the thickness of the fluid and the characteristic length scale in the transverse dir- 



ection. The physics of such processes are reviewed in Oron, Davis & Bankoff, 16 
and more recently in Craster & Matar, [l3]. If one includes surface tension in the 
model, the fluid is governed by a degenerate fourth order parabolic equation. Such 
equations, studied by Bernis & Friedman, [s', Bertozzi & Pugh, [o], and Beretta, 



Bertsch & dal Passo |7 still hold many challenges. See, for example, 12123 24] 
for progress on models which include convection, and Becker & Griin, |5], for a 
thorough review of analytical progress. 

There has been recent interest in equations governing the evolution of a thin, 
viscous film on the outer (or inner) surface of a cylinder of radius R rotating with 
constant angular velocity a;, as in Figure [l] Subject to non-dimensionalization, this 
can be modeled by 

(1) ut + {v^{ueee -^uq - sin i9) -^uju)e = 0, 

Of particular interest is the case n = 3, in which one assumes there is no slip at the 
fiuid-cylinder interface. The solution, li, measures the (dimensionless) thickness of 
the fluid, parameterized by the angular variable, ^ Vt = [— 7r,7r). To derive Q, 
one begins with the Navier- Stokes equations, together with a kinematic bound- 
ary condition for the free surface surface tension. Subject to appropriate physical 
scalings and geometric symmetries, we are left with a dimensionless equation with 
just the uo parameter. We refer the reader to, for instance, ^ and ^ for a full 
derivation. 

Though Buchard, Chugunova, & Stephens, |10 , gave a detailed examination of 
([T]) when = 0, our understanding of the case uj ^ remains incomplete. Some 
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Figure 1. The grey region represents the thin film of liquid on 
the surface of the cylinder, which is rotating with constant angular 
velocity uo. We assume symmetry in the z coordinate, so that u is 
a function of t and only. 



partial results on this case have made assumptions that solutions of ([T]) are strictly 
positive, such as in 22 , where Pukhnachov proved, under a smallness assumption, 
positive solutions of exist and are unique. 

In this work, we further study steady state solutions - non-negative^ time inde- 
pendent solutions u{9) of ([T]), satisfying 

(2) (u'^iueee ^uq- sin i9) + uju)e = 
This can be integrated once, to 

(3) u'^ {ueee uq - sinO) uu = q, 

where the constant of integration, is the dimensionless flux of the fluid. In 
particular, we give an elementary proof that when 7^ 0, both weak and classical 
solutions must be strictly positive; see Section [2j 

For a more extended study of the steady states, it is helpful to explore the prob- 
lem numerically. This can provide information on how, and if, u and its derivatives 
develop singularities as the rotation parameter, cj, goes to zero. For such a task, 
one needs a robust algorithm for computing solutions of (|3|. In Section [sj we for- 
mulate and benchmark an iterative spectral algorithm that is easy to implement 
and rapidly converges. 

Our last result, appearing in Section [4j is motivated by Chugunova, Pugh, & 
Taranets, who proved that for n = 3 there can be no positive steady state 
solutions when the rotation and flux parameters satisfy the inequality 



(4) q > ii 



This was refinement of an earlier bound due to Pukhnachov, 21 . We explore the 
(cj,g) phase space and find evidence that this bound may be sharp. Furthermore, 
it appears that this curve is related to a transition between two different regimes 
of steady state solutions. 
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1.1. Notation. In what follows, we shall use the following notational conventions. 
^per(^) is of continuous 27r-periodic functions on ^7 = [— 7r,7r) with k con- 

tinuous derivatives. H^q^{Q) is the set of square integrable 27r-periodic functions 
on ^7 = [— TT, tt) with k square integrable weak derivatives. For a function ^ in one 
of these spaces, we shall denote its mean by ^. For brevity, we will also sometimes 
write 

to indicate that there exists a constant C > such that 

a < Ch. 



2. PosiTiviTY OF Steady States 

Here, we prove that both classical and weak steady state solutions of ([T]) are 
strictly positive provided 7^ 0. Buchard, Chugunova, & Stephens, [TT], demon- 
strated that for = 0, such solutions need not be positive. 

2.1. Classical Steady State Solutions. As noted, the case n = 3 corresponds to 
a no-slip condition at the surface of the cylinder. With this in mind, it is physically 
obvious that, as long as 7^ 0, any steady state solution of ([T]) with n = 3 will coat 
the entire cylinder - that is, there will be no "dry patches" on the surface of the 
cylinder. The following elementary argument shows that this physical intuition is 
correct, and is in fact valid for a broad range of exponents n. 

Proposition 1 (Positivity of Classical Steady State Solutions). Let u G C^q^{Q) 
he a classical steady-state solution of ([T]) with n> 1. //supp {u) 7^ Vt, then = 0. 

Proof. Assume that u vanishes at ^o- Evaluating ([3| at we have (7 = 0. Con- 
sequently, 

(5) u^{ueee -\- ug — sin 6) -\- uu = 0. 
On the support of li, 

(6) u^'-^iueee ue - sinO) ^ uj = 0. 

Since ueee is bounded as we send 6 towards the root, u = 0. □ 

Though a classical solution of ([T]) needs to be C^, the above proof succeeds 
without change if u is merely with bounded third derivative. 

2.2. Weak Steady State Solutions. Positivity can also be proven under weaker 
assumptions on u. We call a non- negative function u G i^pgr(^) a weak solution 
of([2|,iffor an^GC^r(^). 

(7) / {u''(l)0)e{uee + ^ - sin((9)) - uju(l)e dO = 0. 

J — TT 

We first establish some pointwise properties of a solution in the neighborhood 
of a hypothetical touchdown point, ^o- 
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Lemma 1. Let u G H'^q^{Q) he non-negative. If there exists Oq such that u{Oq) = 0, 
then 

(8a) u{0o) = 

(8b) u{0) = \u{0)-u{0o)\<\0-0of^ 

(8c) \u\0)\ = \u\0)-u\0o)\<\0-0o\'^'. 

Proof. By virtue of a Sobolev embedding theorem in dimension one, ^ C^'^/^, 
|2l[T5]. Therefore, both u and u' are Holder continuous with exponent 1/2; 

(9a) \u{0)-u{9o)\<\9-9o\'/^ 
(9b) |„'(^)_„'(^„)|<|^_^„|i/2. 

Because u is and non- negative, we must have that u'{Oq) = too, otherwise 
there would be a zero crossing. 

Refining (|9|, we can apply the mean value theorem and the Holder continuity of 

to get 

\u{0) - u{0o)\ = \u'{e.)\\e - 0o\ = \u'{0,) - u'{0o)\\0 - 0o\ <\0- 0of^ 
where 6^ lies between 6 and 6. □ 

To prove the weak form of Proposition [l] we begin by deriving a weak analog of 
@: 

Lemma 2. If u is a weak solution in the sense of then there is a constant q, 
depending only on u, such that 

/TV p7T 
{u''(l))0{uoe + ^ - sin((9)) - uju(t)dO = q (t)dO 
-IT J —IT 

for all^eC^,,{n). 

Proof First, observe that for any (j) G C^j,{Q), <P — 4^ has mean zero and therefore 
has an antiderivative in C^^{Q). Given a test function (/), let ip solve cpo = (j) — 4>- 
Plugging in in ([7|) in place of (/), we have 



(11) / {u^(j))e{uee^u-s\ii{6))-uju(j)d6 = (j) I {u'^)e{uee ^u-sm{6)) -uoudO. 

J —IT J —TT 

We may now take 

q = 27T {u^)o{uoo -\- u — s'm{0)) — ujudO. 

J —TT 

□ 

We next derive a weak analog of equation ([5| in the event that u has a root. 

Lemma 3. For n > 1, let u he a weak solution satisfying (10) for all (j) G C^^{Q). 
If u{0o) = for some Oq G ft, then q = 0. 

Proof Given £ > 0, let G C^^{Q) be a compactly supported test function satis- 
fying: 



(12) 0<^<1, supp^ C (-1,+1), / (j)dO = l. 



,+1), [\ 

J —TT 
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Let 

(13) (/>^(^)=£-V((^-^o)/^) 

Then supp (j)^^Is = [Oq — s^Oq ^ e]. Substituting (j)^ into (10), we compute 
1^1 < / u(l)'dO^ [ \nu''-^U0(l)' ^u''(l)0\\uee^u- sin 0\dO 

Jle J Is 

(^^) < \uj\\\u\\l^^i^^ + ||n^i^-^^X0(/)^ + ^i^(/)^|U2(^^)||^i^^ + ^i-sin^|U2(j^) 

< \uj\\\u\\l^^j^^ + CV2~e {\\nu^-^ue^'\\L-^i^) + \Km\L-{h)) 

The constant C depends on the norms of u and its second derivative, which are 
both finite by assumption. 

Using estimates (|8|, we have that within i^, 

(15a) u^4>l<e'^-^ 

(15b) nu"-iue0" <e^-^ 



Substituting this pointwise analysis into (14), 

(16) \q\<ei ^e^^^-^^/\ 

Since this holds for ah e > 0, we conclude that = for n > 1. 



□ 



Using Lemmas |2] and [Sj we now know that a weak solution with a touchdown 
satisfies 

(17) / {(j)u'^)e{uee -^u- sin((9)) - ijou(j)dO = 0. 

J —TV 

Contrast this with its classical counterpart, equation ([5|. 

We now proceed to our main result in the weak case, under the additional re- 
striction that n > 2: 

Proposition 2 (Positivity of Weak Solutions). Let u be a weak solution satisfying 
^ with n>2 for all (p G C^r(^)- If S = {0 : u{0) >0}^n, then cj = 0. 

Proof The proof is by contradiction. Assuming that 5 / 1^, let 6>o be a point on 
the boundary of S. Given e > sufficiently small, take y ^ ^ such that the interval 
h/2 = ~ §7 ^ + |] is contained in S and \0o — y\ = e. Now, let be a compactly 



supported test function as in (12), and let 

r{0) = 2e-'^{2{0-y)/e). 

This function then has support in 1^/2 as in Figure |2] Finally, define (j)^ = tp^ /u^^ 
which is as smooth as ii, by the conditions on e and y. 



Plugging (j)^ into ( |17[ ), we see: 
1^ 

(18) 



u'-''tlj'dO< / ItplWuee ^u- cos 0\dO 

/2 

< We\\LHi,/^)\\ue9 ^ u - cosO\\l2^i^^^). 



Since u ^ e'^ on the interval and n > 1, we can bound the integral on the left from 
below: 

(19) / u^-''ilj'dO>ei^^-''l 

Jle/2 
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By the construction of V^^, ||'^^||l2(/^/2) ^ ^ '^^^ - Therefore, 
(20) |c^| <<r^^"-'^||w^*^* + w-cos^|U.(,^^,) 

For n > 2, we immediately observe that since £ > is arbitrary, = 0. For n = 2, 
we know that since [ueef' and v? are measurable, 

+ ii-sin6>||i^2(j^/2) 

as £ ^ 0, which implies = in this case. 




Figure 2. In the proof of proposition 2, we choose a smooth, 
compactly supported function i\) with support contained in the 
support of li, but close to a zero of u. 



□ 



3. A Computational Method for Steady State Solutions 

While we have given a simple proof that when ^ 0, steady state solutions must 
be positive, we have not established their existence. For that, we rely on 22 , where 
it was shown that in the case n = 3 "small solutions" existed and were unique. Let 
us motivate this with some heuristic asymptotics. First, assume that < g <C cj, 
i.e. 

(21) <7 = ecj, < e <Cl. 

Subject to this assumption, we may then perform a series expansion in e, 

(22) u = u^'^ + eu^'^ + e^u^^^ + e\^^^ + . . . 
Substituting into 

(23) u^{ueee ue — sin 0) -\- uju = euj, 
we then match orders of e. At order e^, 

(24) {u^^^fid^u^^^^deu^^^ - smO)^uju^^^ = 0. 
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As we have no way to solve this, we shah assume the order O(e^) term is trivial, 
= 0; consequently, u is O(e^). Proceeding with this assumption, the next few 
orders are: 

,(1) 



(25a) 
(25b) 
(25c) 



O(e^) 
O(e^) 
0(e3) 



uu^'^'^ = 0, 

,(1)^3 «i 



,(3) 



0. 



(26) 



The leading order approximation of the steady state is thus 

^ , ^ (^\'^ ■ a 
u ^ \ — sm6'. 

UJ UJ \UJ J 

This corresponds to the "strongly supercritical regime" discussed by Benilov et 
al.^ f3^.^ • A similar analysis applies to the case n = 2, for which an approximate 
steady state solution is 

(27) u^^^-m^ ^inO. 

UJ UJ \UJ J 



This regime was studied in 22 , where the expansion was justified by a contraction 
mapping argument. This motivates the numerical implementation of this approach 
as a tool for computing such steady states. 

3.1. Description of the Algorithm. In what follows, we focus on the case n = 3, 
although it can be easily adapted to other cases. We seek to rewrite ^ in the form 

u = L-^f{u)+g{e) 

where L is a linear operator, / contains the nonlinear terms and ^ is a a driving 
term. This will motivate the iteration scheme 

To begin, we divide ([3| by i/^, to obtain 

{dg + de)u = sin^ + qu~^ — uju~'^ . 

As dg + de has a non trivial kernel, the lefthand side cannot be inverted without 
careful projection. Pukhnachov resolved this problem by adding and subtracting 
appropriate terms to the equation. First, we introduce the variable 

(1 

V = u . 

UJ 



Then 



(28) 



Lv ^ {dl 



UJV 



sm( 



UJ 2 3^ + 3gcji' + 'U- 



2^,2 



■Fiv) 



The right hand side of (28) contains a driving term, sin^, and an essentially non- 
linear term of order 0{v^ Rearranging (28), 

(29) V 



We now have our iteration scheme 



(30) 
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In Fourier space, this is 




(a) Solutions with fixed ou = 0.1 (b) Solutions with fixed q = 0.01 

Figure 3. Steady state solutions with various uj and q values. 
Computed by the iterated spectral algorithm at 512 grid points to 
a tolerance of 10~^^. 
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Figure 4. Convergence to a steady state, with v^^"^ 
0.0192 and q = 0.01. 



0, UJ 



3.2. Performance of the Spectral Iterative Method. When our algorithm 
successfully converges, it does so quite rapidly, as shown in Figure [Sj where we 
plot the error between iterates v'^^^ and the final iterate. We see that the error 
is proportional to e~"^, although it is clear that the rate, a, varies with q and 
UJ. Indeed, as shown in Figure |6j there are choices of these parameters for which 
our algoirthm diverges; this is consistent with the threshold Q. While this will 
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. uj = 0.03, q = 0.001 
LJ = 0.001, q = 10-5 
. UJ = 0.0044, q = 0.0001 
• UJ = 0.0192, q = 0.001 
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Iterate 



Figure 5. error between approximate solutions and the final 
iterate, at 512 grid points, computed to a tolerance of 10~^^. 



1000 




_4 _3 _2 -1 1 2 3 4 

e 

Figure 6. A typical divergent sequence of iterates, caused by the 
selection of cj = 0.1 and q = 0.1. The iterative spectral method 
was run to only 30 iterations - after 50 iterations, the solution had 
a maximum of 4 x 10^. 



be further discussed in Section |4j we belive this is closely related to the transition 
between the "small amplitude regime", and an essentially nonlinear regime. 

However, when we do observe convergence, there is a very high degree of preci- 
sion, even using a modest number of grid points. Figure [Tj shows how quickly the 
maximum and minimum of the numerical solution converge to one with 10^ grid 
points as the number of grid points increases. Again, solutions with values of q and 
UJ closer to Q converge less rapidly; the solution in black is much closer to this 
threshold than the other three. 

3.3. Relation to Other Numerical Methods. Steady state solutions of ([T]) 
have been computed in several ways. Some groups, such as Ashmore, Hosoi, & 
Stone, [3], have run the time dependent equation to steady state. More recently. 
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- 0.87, q = 0.01 
= 0.1,9 = 1 X 10-^ 
= 0.0192, q^lx 1Q-' 
= 0.0192, q = 0.001 
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J = 0.87, q = 0.01 
^ = 0.1,9 = 1x10-5 
J = 0.0192, 9 = 1 X 10-' 
= 0.0192, q = 0.001 
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(a) The error in the maximum of the solution, (b) The error in the minimum of the solution. 



Figure 7. Convergence of the iterative spectral method, as meas- 
ured by an error in the minimum and maximum, as the grid size 
is increased. 



Benilov, Benilov, & Kopteva computed the solutions directly by a finite difference 
discretization with a Newton algorithm, [gI. 

An advantage of our approach is that it converges quickly and is extremely 
easy to program. An implementation of the algorithm in Matlab is available at 
|http : //hdl . handle . net/ 1807/300 15 , 

Our approach bears some resemblance to Petviashvilli's algorithm and the re- 
lated spectral renormalization procedure, which have become quite popular within 
the nonlinear dispersive wave community for computing solitons, For 
such equations, including Korteweg - de Vries and nonlinear Schrodinger/Gross- 
Pitaevskii, solitary wave solutions solve semilinear elliptic equations of the form 



(32) 



V'Q + AQ = /(Q), A>0, 



For power nonlinear it ies, f{Q) = Q^, Petviashvilli's algorithm finds a solution by 
seeking a fixed point of an iteratively renormalized nonlinear equation in the Fourier 
domain, 



(33) 



Q(^+i) = A 



"^A+ |fc|2' 



/(A + |fcp)|Q(^')pdk 

—(J) 

//(QO))Q dk 



Q^-^^ is the j-th approximation, and Aj is a renormalization constant at this iterate. 



Our scheme, (31) lacked any such renormalization constant, Aj. This raises the 



question of whether or not such a constant could improve our approach. Our results 
show this is not the case. Introducing the constant 



(34) 



A. 



E 



-ik^ ik ^ 



we can alter our iteration scheme to be 

+ 1) _ A7 / ^^,1 ~ ^k,- 



(35) 



V 



AJ 



2i 



-ik^ + ik ■ 
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where 7 is a constant which we may alter to stabihze or accelerate the scheme. 

An outcome of our experiments is that not only does our method converge with 
7 = 0, there does not appear to be appreciable benefit with any other value; see 
Table [l] For this reason, we do not employ Aj, and only deem our algorithm an 
iterative spectral method. 



Table 1. Number of iterations required for ||iin+i — '^n|| ^ 10~^^i 
for a variety of initial guesses uq. Computed on 2^^ grid points 
with different values of 7. q = 0.001153 and = 0.2 in all cases. 
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69 


Uo a_ - 

UJ 
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2 


2 
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2 


2 




- ^sech(^) 


38 


18 


10 


3 


9 


17 


39 



3.4. Mass Constraints and Non-Uniqueness. Given values of q and a;, our 

algorithm successfully converges to non-trivial solutions. However, for a given value 
of q and cj, there may be more than one solution of To study this phenomenon, 
it is helpful to introduce the mass of a solution 

(36) M = I u{0)dO. 



Benilov, Benilov & Kopteva, |6 tried fixing cj, and numerically explored the rela- 
tionship between the flux and the mass as they varied. They found non- uniqueness 
of the solutions, in the sense that for a certain range of M, there were multiple 
solutions with different values of q. Alternatively, for a certain range of q^ there 
were multiple solutions with different values of M. See Figure 14 of their work for 



a phase diagram, and see our Figure 10 for a similar diagram. 

This non-uniqueness affects our algorithm. For a given (oJ^q)^ the solution that 
the method converges to may not be the solution we desire, as the following example 



shows. Using, AUTO, 14 , if we seek a solution with M = 1 and u = .09, we will 
obtain the solution pictured in Figure [s] (a). This solution has flux parameter 
q = 0.0114039. 

Suppose we now fix = .09, and continue the solution by lowering the mass. 
As seen Figure [s] (b), there will be a second solution, with the same q^ but slightly 
smaller mass. For q = 0.01140392, the second solution is at M = 0.912782 Thus, for 
a given there may be multiple solutions, each with a different mass- which 

solution does the iterative spectral algorithm converge to? As Figure [s] (a) shows, it 
converges to the solution with the smaller mass. The agreement between the AUTO 
solution and our method's solution is quite good; after the AUTO solution is splined 
onto the regular grid of the spectral solution, the pointwise error is O(10~^). 

The inability of our algorithm to recover the M = 1 solution can be traced to 
the development of a linear instability in the iteration scheme. Let v^^"^ = v -\- p^^\ 



where is a solution, and p^-^^ is a perturbation at iterate j. If we linearize (30) 
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9 M 

(a) Two different solutions with uj = .09 and (b) (M, q) pairs of solutions with uu = .09 
q = 0.0114039 

Figure 8. AUTO finds two solutions by continuation methods 
the same u and but with different mass. The Spectral Iterative 
Method converges to the solution with the smaller mass. 



about we have 

(j+i) _ ^-1 + Qq'^cjov + Aquj'^v'^ + cu^v^ 

(37) 



q^ {q + ujv)^ 

If the spectrum of Ly lies strictly inside the unit circle, the perturbation will vanish, 
and the solution v will be stable with respect to the iteration. However, if there 
are points in the spectrum of this operator which lie outside the unit circle, the 
solution will be unstable to the iteration, and we would expect our algorithm to 
diverge. 

Though we shall not attempt to prove properties of the spectrum of L, we can 
discretize our problem and compute eigenvalues of the induced matrix. Computing 
the eigenvalues of a discrete approximation of Ly at the two AUTO solutions, we 
find that the matrix corresponding to the M = 1 solution does, in fact, have an 
eigenvalue outside the unit circle, Aunstabie = 1.017059; see Figure[9| The discretized 
operator is defined as 

(38) Lv [d^^D^ ^i) diag {F'{y)} 

with the differentiation matrix a dense Toeplitz matrix for band limited inter- 



polants, as in 25 . The value of the unstable eigenvalue converges quite quickly as 
we increase the number of grid points, as shown in Table [2] 

The (M, q) phase space of solutions with fixed uo can be more tortuous than that 
shown in Figure [s] (b). Indeed, as we send ^ 0, a singular limit, there may not 
just be multiple solutions of different mass and the same flux, but also multiple 
solutions of different flux and the same mass. Such an example appears in Figure 
10 These cases are more extensively explored in [g] and in Badali et al, [4]. 
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M = 1 




-1L T , 7 L 

-1 -0.5 0.5 1 

M = 0.912782 




Re. A 

Figure 9. The eigenvalues of the discretized Ly for the two solu- 
tions computed using AUTO with uo = .09 and q = .0114039. Note 
that for the M = 1 solution, there is an eigenvalue outside the unit 
circle at A = 1.017059. Computed with 2048 grid points. 

Table 2. Values of the unstable eigenvalue of Lv for the u = .09, 
q = .0114039 and M = 1 solution as a function of the number of 
grid points. 
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M 

Figure 10. The {M^q) phase space of solutions with rotation 
parameter cj = 5 x 10~^. Note that there can be as many as 
six solutions of different mass and the same flux, and as many as 
three solutions of different flux and the same mass. Computed 
using AUTO. 
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UJ 



Figure 11. The (a;, q) phase space of solutions to ([3| at prescribed 
mass values. As the figure shows, even two parameters are insuf- 
ficient to uniquely determine a steady state solution; where the 
curves of constant mass cross, there will be two or more solutions 
with the same uo and q values. All curves of solutions are con- 
strained by Q, as expected. Also note that at fixed q^ for suffi- 
ciently large cj, the relationship is nearly linear. Computed using 
AUTO. 



4. NON-ExiSTENCE RESULTS 



As mentioned in the introduction, in 12 , the authors prove that solutions to (3) 



do not exist for all {uo^q) pairs. Indeed, they show there are no solutions when (4) 
holds. Interestingly, this relation is stated entirely in terms of q and cj, and makes 
no consideration of the mass, M . 

We now explore, numerically, Q. For fixed mass, we compute the associated 
solution at a variety of uo values, determining q. The results are appear in Figure 



11 As [uj^q) cannot uniquely determine we rely on AUTO and its continuation 
algorithms, rather than our own iterative scheme. 

Examining the figure, we remark on the two features found in each fixed mass 
curve. Far to the right of the non-existence line, the uj-q relationship appears to be 
linear. Indeed, for moderate mass values, the relationship appears to be of the form 



(21). It is in this "small amplitude" regime that our iterative spectral algorithm 
succeeds. 
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As uj continues to decrease, the scaling changes to accomodate Q, and appears 
to fohow a curve hke q (x , for a value of f3 > 3/2. For all masses, the deviation 
from Q is less than an order of magnitude. Thus, there is not an obvious small 
parameter from which to formulate a series expansion. We also remark that there 
are many intersections between the curves of constant mass; the lack of uniqueness 
of a solution for a given (cj, q) is quite common. As uj and q approach this essentially 
nonlinear regime, our algorithm becomes limited. As discussed in the preceding, 
section unstable eigenvalues can appear, precluding convergence to solutions with 
a desired mass, and the rate of convergence slows. 

Turning to Figure [T2j we see the tendency of the solutions to form singularities, 
vertical asymptotes, as u tends towards zero for fixed q. These were computed 
using our spectral iterative algorithm; the asymptotes occurred where the algorithm 
failed to converge in a reasonable number of iterations. As these figures show, it 
is possible to get quite close to the theoretical threshold of Q. This suggests that 
the constraint may be sharp, and this warrants further investigation. Even if the 
constant is not correct, the scaling appears to be sharp, and the inequality plays 
an important role in separating two physical regimes. 
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(a) q=lx 10-3. The mass of the final solution (b) g = 1 x lO"'^. The mass of the final solution 
was 0.3951 was 0.1612 



Figure 12. Thickness of solutions with a grid size of 2 as is 
varied. The dashed lines denote the observed threshold at which 
our spectral iterative algorithm fails and the theoretical threshold, 
Q beyond which no solutions should exist. Computed using the 
spectral iterative method. 



5. Discussion 

In this work, we have given a simple proof of the positivity of steady state 
solutions of a thin film equation with rotation, subject to modest assumptions, 
provided n > 2. It may be possible to employ additional regularity properties of the 
steady state solution to lower this threshold; we only relied on Sobolev embeddings 
in our result. 

We have also formulated a very easy to implement algorithm for computing the 
steady state solutions for a given {uj^q) pair. However, in cases where multiple 
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solutions are known to exist, an instability can develop about one of the solutions, 
precluding convergence to this solution. We conjecture that the preference is related 
to the proximity to the non-existence curve, where the parameters have left the 
small amplitude regime. Finding a way to stabilize the spectral algorithm in this 
regime remains an open problem. 

Finally, our computations suggest further exploration of the non-existence curve 
is warranted. Figures [TT] and [12] suggest that the nonexistence relationship may be 
sharp. Furthermore, the phase diagram indicates that this inequality separates two 
physical regimes. In the first, where uj and q are nearly linearly related, the solution 
is nearly constant, and can be approximated by (26). In order to accomodate the 
nonexistence constraint, the phase diagram bends, and we enter a regime where the 
parameters are nonlinearly related, and the steady state solutions cease to be well 
described by unimodal bumps. 
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